% This script rescales and plots impulse response functions for 
% contemporaneous or news TFP shocks.

% Christoph Gortz


close all
clear all

%% Initialisation 
shockPlot = 'z';             % Shock process to be plotted
exo_shockPlot = 'ez8';        % Exogenous shock to be plotted
shockOther = 'v';             % Other Shock process with trend
exo_shockOther = 'ev';        % Other Exogenous shock with trend
vars = [                    % Variables to be plotted
    'y         ';
    'c         ';
    'inve      ';
    'lab       ';
    'labc      ';
    'labi      ';
    ];
nVars = size(vars,1);
IRFtitle = [                % Title of figures corresponding to variables
    'Output                   ';
    'Consumption              ';
    'Total Investment         ';
    'Total Hours              ';
    'C-Sector Hours           ';
    'I-Sector Hours           ';
    ];

        
%this version has same persistence and shock stdev for all models  (all same as second one)      
dataName = ['Output_twosb_v10b_mobileKNoFItfpNoMe_v3         ' % SGU type model - no nominal rigidities and financial frictions
            'Output_twosb_v10b_mobileKNoFItfpNoMeNR_v3       ' % SGU type model - with nominal rigidities and without financial frictions
            'Output_twosb_v12b_B_aggTFP48_rev2noFIdataNoMe_v3' % Baseline for SGU comparison (with nominal rigidities and financial frictions)
            ];
       
figName = ['IRF_SGU_vs_SGUplusNomRigidities_vs_BaselineSGU_' , exo_shockPlot];

figDimensions = [1, 6];     % Dimensions of figure - needs to be consistent
                            % with nVars
irfLength = 20;             % Length of generated IRFs 


%% LOAD INFO AND DETREND
for iData = 1:size(dataName,1)
    eval([ 'load  ' dataName(iData,:) ';' ])
   
% Load variables' irf data into cell
for iVars = 1:nVars
    eval([ 'irfMean(iVars,1) = {oo_.irfs.' strcat(vars(iVars,:)) '_' exo_shockPlot '};' ])
end

for i = 1:numel(irfMean);
    Variables{i,1} = [irfMean{i,1}];
end
% Load shock's irf data
eval([ shockPlot 'Mean = oo_.irfs.' shockPlot '_' exo_shockPlot ';' ])
eval([ shockOther 'Mean = oo_.irfs.' shockOther '_' exo_shockPlot ';' ])

% Transform series of shock z to isolate zd:
vd =[]; vp = 0;
zd =[]; zp = 0;

for m = 1:numel(zMean)
    % Change of v and z in percent
    vp = vMean(1,m)+vp;
    vd = [vd; vp];
    zp = zMean(1,m)+zp;
    zd = [zd; zp];
end    
vd = vd';
zd = zd';

% Descale stationarised variables: K_t = K_t^~*Z_t
yd = Variables{1,1} + (alfac/(1-alfai))*vd + zd;
cd = Variables{2,1} + (alfac/(1-alfai))*vd + zd;
inved = Variables{3,1} + (1/(1-alfai))*vd;
labd = Variables{4,1};
labcd = Variables{5,1};
labid = Variables{6,1};

% Final Detrended Time Series of IRFs
eval([ 'VariablesDetrended_' dataName(iData,:) '= [yd; cd; inved; labd; labcd; labid;];' ]);

end
%% Plot Detrended IRFs

%% Plot IRFs
figure1 = figure('Color',[1 1 1]); % Generates white background around plot

    for iVars = 1:nVars
        varStr = strcat(vars(iVars,:));
        plotVar = strcat(vars(iVars,:));
        irfTitle = strcat(IRFtitle(iVars,:));
        % Subplots
        subplot(figDimensions(1),figDimensions(2),iVars);
%         for iData = 1:size(dataName,1) % Loop over different datasets of model versions
        eval([ 'irfData = VariablesDetrended_' dataName(1,:) '(iVars,1:irfLength);' ]);
        eval([ 'plot(irfData,''-k'',''linewidth'',1.5);' ]);
        hold on
        eval([ 'irfData = VariablesDetrended_' dataName(2,:) '(iVars,1:irfLength);' ]);
        eval([ 'plot(irfData,''--b'',''linewidth'',1.5);' ]);
        hold on
        eval([ 'irfData = VariablesDetrended_' dataName(3,:) '(iVars,1:irfLength);' ]);
        eval([ 'plot(irfData,''-or'',''linewidth'',1.5);' ]);
%         end
        axis tight
        hold off
        eval( 'title( irfTitle );' );
    end
set(gcf, 'Position', get(0,'Screensize')); % Maximise Figure

